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ABSTRACT 

Spherical (nonrotating) accretion flows with small-scale magnetic fields have 
been investigated using three-dimensional, time- dependent MHD simulations. 
These simulations have been designed to model high-resolution (quasi) steady 
accretion flows in a wedge computational domain that represents a small fraction 
of the full spherical domain. Subsonic and supersonic (super-fast-magnetosonic) 
accretion flows have been considered. Two accretion regimes have been studied: 
conservative, or radiatively inefficient, and nonconservative, in which the heat 
released in magnetic reconnections is completely lost. The flows in both regimes 
are turbulent. They show the flattened radial density profiles and reduction of 
the accretion velocities and mass accretion rates in comparison with hydrody- 
namic Bondi flows. In the conservative regime, the turbulence is more intensive 
and supported mostly by thermal convection. In the nonconservative regime, 
the turbulence is less intensive and supported by magnetic buoyancy and vari- 
ous magnetic interactions. We have concluded that steady, supersonic spherical 
accretion cannot be developed in the presence of small-scale magnetic fields. 

Subject headings: accretion, accretion disks — black hole physics — convection 
— MHD — turbulence 



1. Introduction 



Mass accretion onto a central gravitational body plays an important role in the formation 
and evolution of a large variety of astrophysical objects such as planets, stars, galactic nuclei, 
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galaxies, and clusters of galaxies (see Frank, King & Raine 1992). Depending on the amount 
of the angular momentum £ , which the accreting mass carries at the outer boundary radius 
-Rout; accretion flows can take either a disk-like or spherical-like form. Flows with relatively 
high angular momentum, £ ~ ^/GMR out) where G is the gravitational constant and M is the 
central mass, form centrifugally supported accretion disks (e.g., Shakura 1972). Flows with 
low angular momentum, £ <C \/GMR out , can form spherical or quasi-spherical accretion 
flows at radii R > i\jGM in which the centrifugal force is weak and not sufficient to 
balance gravity (Illarionov & Sunyaev 1975). The theoretical study of spherical accretion 
flows is based on an analytical solution discovered several decades ago by Bondi (1952). 
This solution describes the idealized case of isentropic nonmagnetic accretion flows. Since 
that time, the theory of spherical accretion has been significantly developed so that it now 
includes an understanding the role of different physical mechanisms (such as radiative cooling 
and heating, magnetic field dissipation, thermal and radiation transports, etc.) and more 
realistic inner and outer boundary conditions (for reviews, see Frank et al. 1992; Kato, 
Fukue, & Mineshige 1998). More recent studies of spherical accretion considered convection 
(Markovic 1995), accretion onto a magnetic dipole (Toropin et al. 1999), magnetic diffusivity 
due to turbulence (Shadmehri 2004), and vorticity (Krumholz, McKee, & Klein 2005). 

This paper continues the numerical study of radiatively inefficient spherical (nonro- 
tating) accretion flows with magnetic fields conducted by Igumenshchev & Narayan (2002, 
hereafter IN; also see Pen, Matzner, & Wong 2003 for a similar study). IN have demon- 
strated, with help of three-dimensional MHD simulations, that the effects of a magnetic field 
can significantly modify the structure of Bondi-type flows. They argued that even initially 
weak magnetic fields can produce dramatic changes. The main reason for these changes is 
the local, nonuniform release of the thermal energy during the dissipation of tangled mag- 
netic fields in reconnections. This release suppresses the inward motion of mass and results 
in the development of turbulence that is mainly supported by thermal convection. The dis- 
sipation of a magnetic field in reconnections is compensated by the efficient amplification of 
the radial field component in spherical convergent flows (Shvartsman 1971). IN have devel- 
oped a simple analytical theory of spherical convection-dominated accretion flows (spherical 
CDAFs 1 ), which is similar, with regards to the involved basic physics, to the theory of ro- 
tating CDAFs (Narayan, Igumenshchev, & Abramowicz 2000; Quataert & Gruzinov 2000). 
This theory predicts the flattened radial density profile p oc R^ 1 ! 2 in contrast to the steeper 
density profile in asymptotic Bondi flows, p oc R~ 3 / 2 (see the appendix, eq. [A6]). Because 
of this flattened profile, the mass accretion rate in spherical CDAFs is expected to be signifi- 
cantly lower than the Bondi mass accretion rate in flows with the equivalent outer boundary 



Note that IN used instead the name "convection-dominated Bondi flows" which we do not adopt here. 
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where R m is the inner radius of the flows. These properties make the spherical CDAF solution 
to be a prominent candidate for explaining the phenomenon of dim galactic nuclei contained 
supermassive black holes, including our Galactic center Sgr A* (e.g., Melia 1992; Baganoff et 
al. 2001, 2003; Quataert 2002; Ghez et al. 2003; Ho, Terashima, & Okajima 2003; Soria et 
al. 2006), in which accretion disks are typically not observed. Also, spherical CDAFs can be 
employed to explain the problem of missing isolated neutron stars in our Galaxy (Treves & 
Colpi 1991; Blaes & Madau 1993; Turolla et al. 1994; Belloni, Zampieri, & Campana 1997; 
Toropina et al. 2003; Perna et al. 2003) and the observed properties of isolated stellar mass 
black holes (e.g., Fujita et al. 1998; Agol & Kamionkowski 2002). 

The generality of the numerical results of IN and Pen et al. (2003) were limited in the 
important aspect that no steady or quasi-steady accretion flows were obtained. These authors 
studied transient states of the flows, which originated because of specific initial conditions 
and involved moving outward shocks. These results also suffered from insufficient numerical 
resolution, especially in the innermost region of the flows. In addition, the assumption of an 
initial bipolar magnetic field in IN resulted in the domination of a large-scale poloidal field 
and doughnut-like density distribution near the black hole at the later evolution time (note 
that similar field topology was proposed by Bisnovatyi-Kogan & Ruzmaikin 1974). Unless 
such large-scale poloidal field can be naturally present in some objects, other possible field 
configurations, such as a small-scale field with zero net magnetic flux, can likely be natural 
too in spherical accretion flows. The main goal of the present study is the investigation of 
the role of small-scale magnetic fields in such flows. 

In the numerical aspect, we overcome some of the limitations of the IN's approach 
employing a new simulation design. In this design, we assume a permanent injection of 
mass and magnetic field into the computational domain that allows us to obtain steady or 
quasi-steady accretion flows after performing long-time evolution simulations. The numerical 
resolution has been improved by adopting spherical coordinates and conducting simulations 
in the wedge computational domain, which represents a small fraction of the full spheri- 
cal domain (see Fig. 1). The injected field is assumed to have a small-scale component 
in the form of radially extended magnetic loops with a zero net magnetic radial flux (see 
Fig. 2). Modifications of the simulation technique allow us to investigate two limiting en- 
ergetic regimes: conservative (or radiatively inefficient) and non- conservative (in which the 
heat from magnetic reconnections is completely lost). 

The paper is organized as follows. In Section 2 we describe the simulation technique, 
initial and boundary conditions, and algorithm of the mass and magnetic field injection. 
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Section 3 presents numerical results, and in Section 4 we discuss these results and make 
final conclusions. We reproduce some analytic solutions of hydrodynamic accretion flows, 
including Bondi solution, in the appendix. 



We use the numerical method, which is similar to that used by IN and Igumenshchev, 
Narayan, & Abramowicz (2003). The method solves the system of ideal MHD equations (e.g., 
Landau & Lifshitz 1987), which describe the dynamics of nonself-gravitating plasmas in the 
central gravitational field. Originally, the method employed a nonconservative numerical 
scheme that solves the internal energy equation, which includes the reconnection heat term 
Q (for more details see IN); 



where p is the density, e is the specific internal energy, P g is the gas pressure, and v is 
the velocity. Test simulations have shown that MHD solutions obtained using equation 
(2) conserve the total energy quite poorly because of artificial loosing or gaining energy 
in numerical reconnections. In some our test cases, the relative error of the total energy 
conservation was up to 10%. To solve this problem, the method has been modified by adding 
a conservative option. Using this option, the method solves the total energy equation 



instead of equation (2). Here B is the magnetic induction and q is the total energy flux per 
unit square. Note that in finite-difference MHD schemes the magnetic field is reconnected 
on scales in which the minimum is limited by the gridsize. Typically, this gridsize is much 
larger than the physical reconnection scales in the studied problems. Therefore, the finite- 
difference schemes, including our scheme, can not accurately represent all the details in 
the process of reconnection. For our purposes, however, the provided level of accuracy is 
sufficient. This situation is somewhat analogous to the representation of shocks in finite- 
difference hydrodynamic schemes in which the numerical shock thickness is also limited by 
the gridsize. 

In all our simulations, we have used the ideal gas equation of state 



2. Simulation Technique 



P-t + = - p gV • v + Q, 



(2) 




(3) 



P e = (7 - l)pe 



(4) 



and assumed the adiabatic index 7 = 5/3. 
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This method employs the three-dimensional spherical coordinates (R,9,4>). The com- 
putational domain is limited by a narrow four-facet wedge located at the equatorial plane as 
shown in Figure 1. The domain is extended from R [n to R ont in the radial direction and over 
#o and 4>q degrees in the polar and azimuthal directions. We have assumed 6o = <po = 7r/16 
and the number of angular gridpoints N$ x = 30 x 30. These points are uniformly spaced 
in both 6 and <p directions. Gridpoints in the radial direction are spaced logarithmically 
so that the three-dimensional numerical cells at any radius take an approximately cubic 
shape. This provides the direction-independent local spatial resolution in the simulations. 
The number of the radial gridpoints N R = 303, which corresponds to R ou t/Rm = 10. 



2.1. Boundary Conditions 

We have used three different sets of boundary conditions assumed at the azimuthal, 
polar, and radial boundaries, respectively, of the wedge computational domain (see Fig. 1). 
The periodic boundary conditions for both fluid and magnetic field are applied at the az- 
imuthal boundaries. At the polar boundaries, we use the conditions that no streamlines and 
magnetic lines can go inside/outside through these boundaries. This is achieved by apply- 
ing the reflection boundary conditions for fluid and assuming the continuous tangential and 
reflection normal magnetic components across the polar boundaries. 

At the inner R in and outer R out radial boundaries, we apply the absorbing boundary 
conditions for a fluid. This means that the mass can flow freely out of the computational 
domain, but no mass is allowed to return from outside. Conditions for the magnetic field are 
assumed considering "ghost" boundary zones located on the outside of the computational 
domain. We assume that these zones can contain only a radial magnetic field, whose strength 
is determined from the divergence-free constraint V • B = 0. Numerical tests have shown 
that a mass with a frozen-in tangled magnetic field can freely flow out through these radial 
boundaries without the effects of field or mass accumulation. 

Our absorbing boundary conditions at i?i n qualitatively correctly mimic the conditions 
near the black hole horizon where matter is free falling in the strong gravity field and the 
radial component of the magnetic field dominates the tangential component (e.g., Thorne, 
Price, & MacDonald 1986). In the case of accreting stars with rigid surfaces such as white 
dwarfs and neutron stars, which can also have magnetospheres, the inner boundary condi- 
tion will depend on the radiative efficiency of the flows. Radiatively inefficient plasma will 
probably form a slowly accreting (M <C Meondi) subsonic accretion flows, similar to that "a 
tenuous continuation of the star" discussed by Bondi (1952). If plasma efficiently radiates 
its thermal energy near the stars' surfaces (e.g., Shapiro & Salpeter 1975), the considered 
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absorbed inner boundary conditions can be adequate for accretion flows far away from these 
surfaces. 



2.2. Injection of Mass 

To obtain steady, or quasi-steady, accretion flows, we permanently inject mass into the 
computational domain. The main problem here is to minimize the consequences of the 
interaction of outflows, which can originate during simulations, and the injected mass. We 
have employed an injection algorithm, which has some resemblance to the injection algorithm 
used by Igumenshchev & Abramowicz (1999) in simulations of rotating accretion flows. 

We assume that the mass is steadily injected inside a thin spherical layer with the 
radius i? in j, which is located closely to the outer absorbing boundary at R out - This mass 
is distributed uniformly over 9 and <fi directions and has zero velocity. Under the action of 
gravity, the larger fraction of the injected mass forms an accretion flow. The smaller fraction 
of this mass could escape through R OVLt because of a thermal spread and interactions with 
outflows. 

The injected mass is characterized by the specific internal energy e , which deter- 
mines two regimes of hydrodynamic accretion - subsonic and supersonic. The critical value 

(eo)crit ~ e vir where 

e vir = GM/R mi (5) 

separates these regimes: e > (eo) CT it corresponds to subsonic flows and e < (e ) C rit corresponds 
to supersonic flows (see the appendix). 



2.3. Injection of a Magnetic Field 

The injected mass can also carry a magnetic field that is associated with it. The field is 
injected assuming that only the poloidal component Ag of the vector potential A is nonzero 
in the injected mass; the other two components, A R and A^, are set to zero. In this way we 
can produce B R and components, which obey 



Three different magnetic field configurations have been used in simulations: purely 
(unipolar) radial, purely toroidal, and loop configurations. The former two configurations 
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have been chosen to test our numerical method. The purely radial field is generated at the 
beginning of simulations assuming that 

A e = B^^, (7) 

where B m j is the magnetic induction at i?i n j. This field remains unchanged during simulations 
by virtue of the conservation of the radial magnetic flux confined in the wedge computational 
domain. 

Two other field configurations, purely toroidal and loop fields, are formed by perma- 
nently injecting the corresponding field at R m j. This field injection is tightened to the mass 
injection and performed by correcting Ag inside the injection layer at each time step as 
follows: 

A' e = A e + £AA, (8) 

where A A is an increment and £ is a correction factor < £ < 1. The increment A A depends 
on the assumed field configuration. In the case of purely toroidal field, 

AA = B R mi , (9) 



where B = \/8n / y( / y — l)e Ap/(3 and Ap is the density of the injected mass. In the case of 
loop field, 

A.4 = BAj %l£^, (10) 
R m 

where the parameter m determines the number of azimuthal sectors in which the radial 
magnetic field periodically changes direction. A combination of two such sectors forms a 
magnetic loop. In the present simulations, we have assumed two magnetic loops in the 
computational domain as illustrated in Figure 2. 

The correction factor £ in equation (8) is used to maintain the strength of the injected 
field at a given level determined in terms of the plasma (5 = P g /P m , where P m = B 2 /8n is 
the magnetic pressure. We assume that 

</3>-/V 



£ = min 



1 , max I 



(11) 



A> 

where (3q is a parameter, which determines the field strength, and (f3) is the (5 averaged over 
the volume of the injection layer. 



3. Numerical Results 



We initiate simulations from a static, nonmagnetic, very low-density medium that fills 
the entire computational domain. Accretion flows are created by steadily injecting the 
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mass into this medium at R m j. These flows generally go through an initial transient phase 
before relaxing to a steady or quasi-steady state. The transient phase typically takes a few 
tens (~ 10 — 30) of the free-fall time, ~ 0.67-R^/VGM, and includes the formation and 
dissipation of shocks and nonlinear waves. The hydrodynamic flows and MHD flows with the 
purely radial and toroidal fields have finally been relaxed to steady states. The MHD flows 
with the loop field have remained time variable, or quasi-steady, even after the completion 
of the transient phase. 

IN discussed in detail the role of magnetic reconnections in spherical accretion flows 
with f3 ~ 1. This role depends on the amount of reconnection heat contributed to the gas 
internal energy. We have considered two extreme energetic regimes, conservative and non- 
conservative. In the conservative regime, the dissipated magnetic energy is fully transformed 
into the heat. Models in this regime have been calculated by employing the total energy 
equation (3). In the nonconservative regime, the dissipated magnetic energy is totally lost 
and does not produce any heat. In this regime, we employ the internal energy equation 
(2) in which the dissipation term Q has been set to zero. Note that the term descibed the 
artificial heat from shocks has been retained in equation (2). 

The study of these conservative /nonconservative accretion regimes can be astrophysi- 
cally motivated. The conservative regime can corresponds to the very high accretion rate 
flows with M ^> MEdd, which do not radiate much energy due to the large optical depth of 
the flows (e.g., Katz 1977). Here Mgdd — 4:irGM/aTm p c is the Eddington accretion rate, 
ot is the Thomson-scattering cross-section, and m p is the proton mass. Another option for 
the conservative regime is the very low accretion rate flows, M <C MEdd, if one assumes that 
the energy released in magnetic reconnections primary does to the ions, which cannot loose 
this energy efficiently via the usually considered electron-ion Coulomb collisions because of 
the tenuous plasma (e.g., Esin et ai). There is another possibility, however, that a signif- 
icant fraction of the reconnection energy will go to electrons, which can efficiently radiate 
this energy by variety of mechanisms (Bisnovatyi-Kogan & Lovelace 1997, 2000; Quataert & 
Gruzinov 1999). The latter possibility provides the additional cooling mechanism for ions, 
which in some circumstance could be more efficient than the Coulomb collisions, making the 
low, or even very low, accretion rate flows to be nonconservative, or radiatively efficient. 

We shall describe models that are either subsonic or supersonic and differed by the 
strength of the injected fields, amount of the heat release in reconnections, and field con- 
figuration. The model parameters are listed in Table 1. For convenience, we have divided 
all these models into three groups: Bondi-type flows (Models A1-A4), subsonic MHD flows 
(Models Bl and B2), and supersonic MHD flows (Models CI and C2). The structure of the 
Bondi-type flows, i.e. flows without magnetic fields or with the fields of special topology, is 
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similar to the structure of Bondi (1952) solution and well approximated by nonmagnetic ana- 
lytic solutions discussed in the appendix. The structure of the latter two groups of models is 
significantly modified because of magnetic fields. We shall describe these groups separately. 

3.1. Bondi-type flows 

Bondi-type flows (Models A1-A4, see Table 1) have been calculated to verify our numer- 
ical method and obtain reference models. These flows are steady and laminar. Models A1-A3 
represent subsonic flows and Model A4 represents a supersonic flow. Models Al and A4 have 
no magnetic field, and thus describe hydrodynamic flows. Model A2 has the purely toroidal 
but sufficiently subequipartition ((3 ^> 1) magnetic field and Model A3 has the uniformly 
ordered strong (/3 < 1) radial magnetic field. In the latter two cases, the magnetic field does 
not affect the structure and dynamics of the flows and these cases are almost equivalent to 
Model Al. 

We shall describe the evolution and structure of subsonic Bondi-type flows using Model Al 
as a representative example. This model has evolved from the initial state (see §2) to a steady 
state through the transient phase. This phase includes the formation of temporal shocks and 
waves, which many times propagate inward and outward in the radial direction, reflecting 
from the inner and outer boundaries, before complete disappearance after t ~ 20 — 30t s . The 
final flow in Model Al is steady, effectively one-dimensional, and described by the analytic 
solution (eqs. [A4] and [A8]). This solution suggests that Model Al represents a part of the 
flow located deeply inside the Bondi radius, R out = 0.0371/?^, and, therefore, this model 
can be closely approximated by an asymptotic Bondi solution (A6), which appears in the 
limit R <C Rb- Figure 3 shows selected properties of Model Al (short-dashed lines) as func- 
tions of the radius. These properties include the distribution of density p (left top panel); 
gas pressure P g (left middle panel); relative temperature T/T vir = e/(GM/R) (left bottom 
panel); Mach number M. — v/c s (right top panel), where c s = \J^jPj p is the sound speed; 
relative radial velocity v/v s (right middle panel); and relative accretion rate (right bottom 
panel). In the latter case, the accretion rate is normalized to the mass injection rate (see 
§2). Because the flow is steady, the accretion rate is independent of the radius. Note that 
about 90% of the injected mass forms inflow. The other 10% escapes through R out because 
of a thermal spread. 

For comparison, Figure 3 also shows a self-similar solution (A9) in which Ai = 1 (long- 
dashed lines). This solution is "boundary free" and has asymptotic Bondi profiles p oc R~ 3 ^ 2 
and P g oc R~ 5 ^ 2 throughout. Note that Model Al demonstrates slightly flatter density and 
pressure profiles, which can be explained by the deviation of R/ Rb from the zerolimit. 
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Model A2 has the weak and dynamically unimportant toroidal magnetic field > 1. 
In all other aspects, Models A2 and Al are similar. Using Model A2, we have tested the 
ability of our numerical method to accurately model the passive transport of magnetic field. 
The magnetic flux conservation requires that the toroidal field near the equatorial part of a 
spherical flow is changed as 

P oc (Rv)- 1 . (12) 

This leads to oc R~ x l 2 and P m oc R" 1 for the self-similar velocity profile v oc R" 1 ! 2 . 
Because Model A2 has the steeper velocity profile, however, the actual profile of P m obtained 
in this model is accordingly flatter and its slope is fully consistent with estimate (12). 

Model A3 has a uniform radial magnetic field. This field is changed with the radius as 
Br oc R~ 2 and therefore P m oc R~ A . The strength of the field has been chosen such that it 
has (3 = 30 at P out and /3 ~ 0.1 at Pq n . The simulations have shown no effects from the flow 
transition from the subequipartition to superequipartition field regions. We shall see later, 
however, that such a transition causes significant changes in the flows with the loop field. 

Note the artificial nature of the field topology assumed in Model A3 that represents a 
small sector of the spherically symmetric monopole field. Such a field can not have magnetic 
reconnections, which actually play an important role in realistic MHD flows. Accretion flows 
in a strong unipolar magnetic field similar to that considered in Model A3, however, can 
occur at magnetic poles of neutron stars and white dwarfs. 

Properties of the supersonic hydrodynamic Model A4 are shown in Figure 4 with short- 
dashed lines. The density and pressure profiles in this model demonstrate nonmonotonic 
behavior at the outer region and approaching the asymptotic Bondi powerlaws p oc P~ 3 / 2 
and P g oc Rr h l 2 at the inner region (see left top and left middle panels in Fig. 4). The 
Mach number and accretion velocity (see right top and right middle panels in Fig. 4) show 
a significant increase inward from the outer boundary. All these properties of Model A4 are 
in good agreement with analytic solution (All), which fits Model A4, assuming P out /P_B = 
0.071 (see the appendix). 



3.2. Subsonic MHD Flows 

Subsonic MHD flows are represented by Models Bl and B2 (see Table 1) and a have 
dynamically important magnetic field (3 ~ 1 — 10. This field is injected into the computational 
domain in the form of magnetic loops stretched in the radial direction (see Fig. 2). The 
strength of the field is determined by the parameter (3q (see Table 1 and eq. [11]). The total 
energy is conserved in Model Bl, whereas the energy released in magnetic reconnections is 
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completely lost in Model B2. 

Simulations of Models Bl and B2 have been initiated from the steady hydrodynamic 
Model Al by gradually increasing the strength of the injected magnetic field to provide a 
smooth transition from nonmagnetic to magnetic flows. Even in this case, however, Mod- 
els Bl and B2 have been settled into their final quasi-steady states after passing through 
the initial transient phases. The final quasi-steady states are turbulent and characterized 
by random variations of all quantities on all spatial scales from the grid size to the size of 
the computational domain. The amplitude of these variations, however, is not large and 
the time-averaged properties of the accretion flows remain unchanged on the large time 
scale. The turbulence in Models Bl and B2 is clearly originated because of the effects of the 
magnetic field. 

Figures 5-9 illustrate the structure of conservative Model Bl, showing the snapshots 
of two-dimensional distributions of selected quantities in the equatorial cross section of the 
three-dimensional computational domain. Figure 5 represents the density distribution, which 
shows clearly recognizable small-scale (with respect to the scale R) density variations of the 
relative amplitude ~ 2. The lower and higher density regions typically take the form of 
filaments, which are predominantly extended in the radial direction. These radially extended 
structures can also be seen in the velocity streamlines in Figure 6. The streamlines form 
a complicated flow pattern consisting of the radially extended narrow inflowing/outflowing 
streams and small-scale vortices. 

The time- dependent flow pattern in Model Bl demonstrates the randomly repeating 
events of the interchange instability. In these events, colder, denser matter is accumulated 
above (i.e., far from the center) the region with hotter, lower density matter. With time, this 
denser matter begins to move down through the low-density region, forming a characteristic 
radially inflowing dense stream. Such a stream typically forms a "mushroom" at its head 
and propagates about half a radius inward from the radius of origin. At the same time, 
the stream carries a frozen-in magnetic loop whose field strength is amplified because of a 
radial convergence of the stream. The example snapshot of magnetic lines in Figure 7 shows 
many such magnetic loops resulting from the interchange instability. During the subsequent 
evolution, the inflowing radial streams are fragmented into small pieces. This fragmentation 
is typically triggered by reconnections of the oppositely directed magnetic lines confined in 
the streams. The reconnections locally release energy and produce a hot low-density matter. 
This matter has a positive buoyancy and forms narrow low-density outflows. 

Figure 8 shows the snapshot of the plasma (5. The spatial distribution of (3 and, ac- 
cordingly, magnetic field is highly nonuniform. The small-scale, high-/? regions, (3 > 100, 
cover the entire computational domain and correspond to the weak field regions, which are 
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typically associated with reconnection regions. The large number of these regions clearly 
indicates the high efficiency of the reconnection dissipation. The low-/? regions, /3 < 1, 
or, accordingly, the strong field regions are typically elongated in the radial direction and 
associated with the inflowing streams. The field in such regions becomes especially strong, 
/3 <C 1, in the innermost part of the flow, R < 2R m , where the high- velocity inflowing streams 
dominate in the flow pattern. 

The chaotic inward and outward turbulent motions in Model Bl are only partially 
supported by magnetic interactions, whereas the more important support is provided by 
thermal convection. As discussed earlier, magnetic reconnections are very efficient in the 
considered model and produce hot, low-density, narrow outflows. Sometimes several such 
outflows can coalesce and form a large convective bubble, which moves outward easier and 
faster. The representative example of such a bubble can be seen in Figure 5 as the low-density 
region near the outer boundary and, in Figure 9, as the corresponding increase of the specific 
entropy. Later in time, this particular bubble has escaped through R ont - Figure 9 also shows 
many other local regions of high specific entropy. These regions typically coincide with the 
lower density regions in Figure 5, have positive velocities, and can therefore be identified as 
convective bubbles. However, not all such bubbles will be able to escape through R ont - A 
significant part of these bubbles will be pulled inward, mixed with the inflowing cold matter, 
and absorbed at R m . 

To make the argument concerning the development of convection more quantitative, 
we have calculated a one-dimensional time- and angle-averaged distribution of the specific 
entropy in Model Bl (see the description of the averaged procedure below). This distribution 
has the negative slope and therefore satisfies the Schwarzschild criterion for convection. 

The nonconservative Model B2 is similar in many respects to the conservative Model Bl. 
Model B2 also forms a turbulent quasi-steady flow pattern. Like in Model Bl, the devel- 
oped nonuniformities in the density frequently take the form of narrow filaments, which are 
extended in the radial direction. However, the amplitude and intensity of turbulent fluctu- 
ations in Model B2 are significantly reduced in comparison with Model Bl. As the result, 
the rate of reconnection dissipation is smaller in Model B2 and the average magnetic field is 
stronger. The turbulent motions in Model B2 involve different MHD effects and are mainly 
supported by magnetic reconnections, magnetic buoyancy, and interchange instability. Con- 
vection motions are not developed in this model because of the absence of reconnection heat. 
The latter explains the less efficient turbulence observed in this model. 

Our MHD models are time variable and, therefore, it is more practical to study the 
structure of these models by performing time and space averaging. In this way, we have 
constructed one-dimensional distributions of different quantities, averaging them over the 
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9 and directions. We have also averaged the obtained spatially averaged distributions 
over time, assuming the time- averaged interval r ~ 3tg. Figure 3 shows the averaged radial 
profiles of selected quantities in Models Bl and B2 (solid and dotted lines, respectively). A 
comparison of these profiles shows certain differences in all quantities except in the accretion 
rates (see lower right panel in Fig. 3). These rates are close to each other and almost 
independent of the radius. The latter confirms that the considered flows are in quasi-steady 
states. 

Model B2 shows a little steeper average profiles for p and P g and a flatter profile of v in 
the inner and middle radial ranges, R < 5i?i n , in comparison with Model Bl. The profiles of 
P m show almost the same slopes in both models, but the magnitude of P m is by the factor 
of ~ 2 larger in Model B2 than in Model Bl. As discussed earlier, this difference in P m can 
be explained by the less efficient dissipation of magnetic field in Model B2. 

The presence and absence of the reconnection heat resulted in different temperatures in 
Models Bl and B2 (see left lower panel in Fig. 3). In agreement with a naive expectation, 
the conservative Model Bl has the larger temperature. At the same time, this model has 
a larger magnetic Mach number, A4 m = v/ ^ c 2 s + c A , where ca = a/ B 2 /Anp is the Alfvenic 
speed, and larger v in the inner region (see right top and right middle panels in Fig. 3). The 
latter two facts can not be simply understood using the analogy with Bondi-type flows in 
which lower Ai and v correspond to hotter flows. It seems that the stronger magnetic field 
is responsible for the relative reduction of these quantities in Model B2. 

The averaged radial structures of MHD Models Bl and B2 show noticeable differences 
from the radial structures of hydrodynamic Model Al and self-similar solution (A9) (see 
Fig. 3). The MHD flows have the flattened density profiles. In terms of the power-law 
distribution p oc R~° ', Models Bl and B2 have a ps 0.7, whereas Model Al and solution (A9) 
have respectively a ~ 1.3 and o = 1.5 (see left upper panel in Fig. 3). Similarly, Models Bl 
and B2 have the flattened profiles of P g oc i?"* 7 in which r] ~ 1.5, whereas Model Al and 
solution (A9) have rj 2.2 and r\ = 2.5, respectively (see left middle panel in Fig. 3). 
Other differences between hydrodynamic and MHD flows include the larger temperatures 
and reduced accretion velocities in the latter flows (see left lower and right middle panels in 
Fig. 3). The accretion rates in Models Bl and B2 show an ~ 10% reduction with respect to 
Model Al (see right lower panel in Fig. 3). This is a clear indication that turbulent MHD 
accretion flows have reduced accretion rates. 
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3.3. Supersonic MHD Flows 

Supersonic MHD flows are represented by Models CI and C2 (see Table 1). These 
models have been initiated from the supersonic hydrodynamic Model A4 by injecting a loop 
magnetic field at i?i n j. Similar to the subsonic MHD models (see §3.2), Models CI and C2 
have approached quasi-steady states after passing through the transient phases. However, 
the transient phases in these models show some differences: the developed temporal shocks 
and waves are downshifted by the supersonic inflow and never reached the outer boundary 
R OU f in quasi-steady states, each of these models obtains a new feature — a quasi-steady 
shock. This shock divides the flows into two regions: the outer and inner, which have super- 
fast and sub-fast-magnetosonic accretion velocities, respectively. In the outer region, the 
flows are radial and laminar. The magnetic field is dominated by the radial component; this 
is quickly increased inward as B oc R~ 2 . No significant field dissipation and reconnections 
have been observed in this region. In the inner region, on the other hand, reconnections 
become important and the flow is turbulent. The quasi-steady shocks in both Models CI 
and C2 can be classified as Alfvenic shocks; they change the orientation of magnetic lines, 
but leave the field strength almost unchanged. Right after the completion of the transient 
phases, these shocks have been located at R ~ 2R in in both Models CI and C2. 

During the subsequent evolution, Models CI and C2 have been slowly evolved on the 
time scales ^> t s . As a result of this evolution, the inner regions in these models have been 
expanded, and, consequently, the Alfvenic shocks have slowly moved outward. The shock 
moves relatively faster in the case of the conservative Model CI and slower in the case of the 
nonconservative Model C2. Specifically, during the evolution time t ~ AOtg counted from the 
end of the transient phases, the shock has propagated a distance of ~ 6 Ri n in Model CI and 
the distance ~ i?i n in Model C2. It is worth noting that these slowly moving Alfvenic shocks 
do not have analogy in hydrodynamic accretion flows. In the latter flows, any developed 
shocks are nonstationary and move either inward or outward, depending on the assumed 
conditions, on the timescale of ~ t ff . 2 

Figure 10 illustrates the flow pattern in the conservative Model CI, showing the velocity 
streamlines projected onto the equatorial plane. The Alfvenic shock is located at R 7i? in 
in the shown moment. Two flow regions can be clearly distinguished: the pre-shock laminar 
outer and post-shock turbulent inner regions. Figure 11 shows magnetic lines that correspond 
to the flow pattern shown in Figure 10. In the pre-shock region, magnetic lines are purely 
radial, but oppositely directed in different sectors (see Fig. 2). These lines do not experience 
reconnections because of the purely radial flow pattern. In the post-shock region, magnetic 



2 We do not consider the radiative shocks here, which can be (quasi) steady. 
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lines are tangled because of the turbulent flow pattern and frequently reconnected. 

The difference in Alfvenic shock motions in Models CI and C2 suggests that there are two 
different mechanisms that cause these motions. Model CI, which includes the reconnection 
heat, develops a convection in the post-shock region, similar to the convection in Model Bl. 
This convection transports the heat outward, causing the relatively fast expansion of the 
inner region and faster motion of the shock. This mechanism, however, is not suitable for 
Model C2, which does not include the reconnection heat and does not develop convection. 
Instead, the shock motion in Model C2 can be explained by magnetic buoyancy, which acts 
in the post-shock region in a manner similar to the action of magnetic buoyancy in Model B2. 
As a result, the magnetic field and energy are transported outward, forcing the inner region 
to expand and move the shock outward. The slower expansion of the post-shock in Model B2 
is explained by the less-efficient energy transport provided by magnetic buoyancy. 

The radial distributions of selected quantities in Models CI and C2 are shown in Figure 4 
by solid and dotted lines, respectively. These distributions have been obtained employing 
the averaging techniques described in §3.3 in which the time-averaged interval is chosen to be 
smaller, r ~ O.ltg, to avoid "washing out" the moving shocks. One can clearly see that the 
structure of Models CI and C2 changes quite sharply at the Alfvenic shocks, whose locations 
at the shown moments are R ~ 7.5 P; n and 3 i?i n , respectively. 

The outer regions in Models CI and C2 are almost identical to each other and similar to 
the outer part of Model A4 (the dashed lines in Fig. 4). This can be seen in the distributions 
of all quantities except for the magnetic Mach number M m (see right upper panel in Fig. 4). 
Models CI and C2 include strong magnetic fields that significantly reduce M. ni (but still 
Ai m > 1 before the shocks), whereas the accretion velocities remain almost the same as in 
hydrodynamic Model A4 (see right middle panel in Fig. 4). Note the interesting behavior 
of P g and P m in these models (see left middle panel in Fig. 4). The values of P g before the 
shocks in the MHD models closely follow the corresponding value in Model A4. The values 
of P m , which are started from the subequipartition level at P; n j, (3 = 10, quickly exceed 
the equipartition level, increasing inward as P m oc R~ 4 . In the case of Model C2, (3 0.1 
just before the shock and the energy of the flow is dominated by the kinetic and magnetic 
energies, which are of approximately equal magnitudes. 

The inner turbulent regions in Models CI and C2 are sub-fast-magnetosonic, M m < 1, 
and have Alfvenic or moderately super- Alfvenic accretion velocities, v > c a . These regions 
are more dense, more hot, and have lower accretion velocities than the corresponding part of 
Model A4 (see Fig. 4). Magnetic reconnections are important here; they support turbulence 
and provide dissipation of the magnetic fields. This dissipation noticeably reduces the slope 
of P m with respect to the corresponding slope in the outer region of the flow where the 
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dissipation is negligibly small. In the case of Model C2, the power-law index r] for the 
magnetic pressure distribution P m oc R~ v is changed from rj = 4 in the outer part to 77 ~ 3 
in the inner part of the flow (see left middle panel in Fig. 4). 

4. Discussion and Conclusions 

We have numerically investigated quasi-steady MHD spherical accretion flows with im- 
posed small-scale magnetic fields. We have confirmed the previous theoretical expectation 
and numerical results such that the flows are turbulent and have radial structures different 
from the Bondi-type accretion flows (e.g., Shvartsman 1971; IN). Turbulence in the MHD 
flows is developed and supported by the interchange instability, thermal convection, and 
various magnetic interactions, including magnetic reconnections and buoyancy. The highly 
nonuniform release of energy in reconnections, which is, in fact, the release of the gravi- 
tational energy converted and stored in the form of magnetic field, makes these flows so 
special and different from the laminar and stable Bondi-type flows (however, see Kovalenko 
& Eremin 1998). We have found that magnetic buoyancy, in addition to thermal convec- 
tion, can play an important role in the modification of the flow structure, especially in the 
case of nonconservative (or radiatively efficient) flows in which the convection could not be 
developed. 

The most important consequence of turbulence in our MHD models, both conservative 
and nonconservative, is the modification of the radial flow structure. Figure 3 demonstrates 
this modification. Turbulent subsonic MHD Models Bl and B2 have the flattened time- 
averaged density profiles, higher temperatures, and lower accretion velocities in comparison 
with their hydrodynamic counterpart Model Al. These properties make Models Bl and 
B2 more like spherical CDAFs (see IN) than Bondi accretion flows. The theory of spherical 
CDAFs predicts the flattened power-law density profile p oc R~~ a in which a = 0.5. Models Bl 
and B2 have profiles close to this, but more steep, a 0.7. These steeper profiles can be 
explained by the influence of the inner boundary condition in our numerical models, whereas 
the analytic CDAF solution is "boundary free." Abramowicz et al. (2002) argued that the 
proximity of the black hole absorbing boundary makes the inner regions of CDAFs to be 
advection dominated. In particular, they demonstrated that viscous rotating CDAFs become 
advection dominated inside ~ 50 to 100 gravitational radii. Our numerical models have a 
rather small radial range, R out /Ri n = 10, and, therefore, the effects of the inner boundary 
can be significant. 

Turbulent Models Bl and B2 have shown the reduction of the accretion rates in com- 
parison with laminar Model Al (see right lower panel in Fig. 3). This result qualitatively 
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confirms the predicted reduction of the accretion rate in spherical CDAFs (see eq. [1]). 
However, the actual reduction of the accretion rates, which is about 10%, shows a poor 
quantitative agreement with estimate (1). We explain this poor agreement by the limited 
radial range in our models, whereas estimate (1) was obtained in the limit R out ^> R[ n . 

The conservative and nonconservative subsonic MHD models (Models Bl and B2, re- 
spectively) are turbulent, however, the turbulence properties in these models are different. 
In the conservative model, the turbulent motions are more intensive and mainly supported 
by thermal convection, which makes this model similar to spherical CDAFs. The noncon- 
servative model has less intensive turbulence, which is supported through various magnetic 
interactions in which magnetic buoyancy seems to be dominated. In some respect, the mag- 
netic buoyancy acts similar to the thermal convection; it also transports (magnetic) energy 
outward. This transport can explain the flattened density profile in the nonconservative 
model, like the convection in spherical CDAFs (see IN). However, the energy transport pro- 
vided by magnetic buoyancy is less efficient and much weaker than the convection transport. 

Our supersonic MHD accretion flows (Models CI and C2) have been initiated by in- 
jecting a low entropy matter into the computational domain. These flows form quasi-steady 
Alfvenic shocks that separate the laminar super-fast-magnetosonic outer inflows and the 
post-shock turbulent nearly Alfvenic inner inflows. We have found that the averaged flow 
pattern in these models is not steady on the large time scale ^> tg. The post-shock regions 
gradually expand, forcing the Alfvenic shocks to move outward. These regions are quite 
similar to the subsonic MHD flows in Models Bl and B2, and, therefore, the gradual expan- 
sion of these regions can be explained by the outward energy flux provided by convection 
and/or magnetic buoyancy. We have found that the Alfvenic shock moves significantly faster 
in the case of Model CI, which develops a convection in the post-shock region. The latter 
is consistent with our observation that the more intensive turbulence and, respectively, the 
larger outward energy flux happens in convective Model Bl, rather than in nonconvective 
Model B2. 

Based on our study of the supersonic MHD models, we conclude that stationary su- 
personic (or super-fast-magnetosonic) accretion flows cannot be realized in the presence of 
small-scale magnetic fields. This conclusion is equally applied to radiatively efficient and 
inefficient flows. Such supersonic flows will unavoidably create shocks at the equipartition 
radius and these shocks will be moved outward because of the action of convection and mag- 
netic buoyancy, which are developed in turbulent post-shock regions of these flows. These 
post-shock regions will continue to expand and, on the large time scale, fill the entire accre- 
tion domain, causing the flows to be subsonic everywhere. 

The wedge computational domain of our models has the limited opening angle (see 
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Fig. 1), and we use the specific boundary conditions in the angular directions (see Sec- 
tion 2.1) to minimize the influence of this limited-size domain. Definitely, the used boundary 
conditions could have some effects on the properties of the simulated MHD turbulence, for 
example, limiting the spatial scales and effecting the motions in the vicinity of the polar 
sliding boundaries. Other consequence of the employing this domain is the inability to sim- 
ulate large scale magnetic structures, which can be developed during the reverse cascade of 
energy from small to large spatial scales in a MHD turbulence. This issue of the limited 
computational domain should be addressed in future works. 

The accretion of mass in our MHD models is accompanied by the reconnection dissipa- 
tion of the magnetic field in turbulent flows. The rate of dissipation is consistently regulated 
through feedback mechanisms that also regulate the intensity of the turbulence. To illus- 
trate the dissipation process quantitatively, we consider the nonconservative MHD models 
in which the change of the time-averaged radial flux of the total energy 



R 



R 2 q R dtt (13) 

n 



directly corresponds to the amount of dissipated magnetic energy. The integral in equation 
(13) is taken over the solid angle of the computational domain Q and the radial flux per 
unit square (see eq. [3]) is 

(v 2 P B 2 GM\ B R . „, , 

te = ^_ + e+ _ + ___j__£ (v .B). (14) 

The solid lines in Figure 12 represent the radial dependence of the total energy fluxes obtained 
in Models B2 and C2. These fluxes have been normalized to the flux GMM in /R in . The 
dashed lines in Figure 10 correspond to the flux conservation Fr = const. The difference 
between the dashed and solid lines at a given R represents the amount of the energy dissipated 
in the radial range from R ont to R. Note that the energy dissipates in the whole volume in 
Model B2, whereas it dissipates only in the post-shock region and no energy dissipates in 
the pre-shock region in Model C2. The total amount of reconnection losses given in units of 
gravitational energy at _R in is about 6.5% in Model B2 and about 5% in Model C2. 
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A. Some Analytic Solutions of Spherical Hydrodynamic Accretion Flows 

In this paper, we have compared our numerical models and analytic solutions of spherical 
hydrodynamic accretion flows. We reproduce some relevant formulas for these solutions 
below. 

Spherical, stationary, and isentropic hydrodynamic accretion flows are described by the 
continuity equation 

M = 4irR 2 pv, (Al) 

Bernoulli's equation 

v 2 7 P GM 

Y + ~ -R- = instant, (A2) 

and the politropic equation of state 

P = Kp\ (A3) 

where P is the gas pressure, K is the politropic constant, 7 is the politropic index, and M is 
the accretion rate. Bondi (1952) solved these equations assuming that the flow is at rest and 
of the uniform density p,^, pressure P^, and sound speed Cqo = a/7-Poo/Poo a ^ infinity The 
Bondi solution depends on the spatial scale R B = 2GM/c 2 O0 , called the Bondi or accretion 
radius, and can be expressed in the following implicit form that gives the radial dependence 
of p: 

Rb 



R J \ p J 7-I 



> 7-1 

Po 



R 0, (A4) 



where A = M/(47ri?|p 00 c 00 ) is the dimensionless accretion rate defined by the expression 

57 — 3 37 — 5 

/ ] \ 2( 7 -l) / 5 _ 2( 7 -l) 

A =U) {-rj ■ < A5 > 

In the small-radii limit, R R B , the Bondi solution (A4) is represented by the asymptotic 
power-law distributions 

v oc R- 1/2 , 

p oc iT 3/2 , (A6) 
P oc R- 5 ' 2 . 

In the case ofl<7<5/3, the Bondi solution (A4) describes transonic flows in which the 
Mach number 

/ p \ 2 / n \ (7+l)/2 

M = \(*£) (A7) 



R J V P 
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is zero at infinity and monotonically increases inward, approaching the infinite value at 
R = 0. 

In the subsequent discussion, we shall consider only the case 7 = 5/3, which is known to 
be the special case of the Bondi solution (A4). In this case, the flow is subsonic everywhere 
and M. — 1 only in the limit of R — > 0. However, such a "boundary free" solution is not 
practical for numerical applications and we have modified the Bondi solution to include 
the finite inner boundary radius R m > in which we require A4 — 1. The modified Bondi 
solution is described by same equation (A4) and the same equation (A7) represents the Mach 
number, but a new A is defined as follows: 

1 + 35 x 1 



where 5 = R in / Rb is a free parameter; < 5 < 1. This modified solution approximates well 
our subsonic Model Al, assuming 5 = 0.00371 (short-dashed lines in Fig. 3). 

In the case of 7 = 5/3, equations (Al)-(A3) allow a self-similar solution, which requires 
the special outer boundary condition = or, equivalently, to set the constant on the 
right hand side of equation (A2) to zero. This solution reads 



v = a^GM/R, 

p = M/(47raVGMR 3/2 ), (A9) 
P = l( 1 _^\ GM P 



5 V 2 ) R ' 

where a is a free parameter; < a < Solution (A9) can be subsonic or supersonic 
depending on a. This solution is characterized by the constant Mach number 



M = ay/l.5/(l -a 2 /2), (AlO) 

and it is supersonic at 1/ \[2 < a < y/2. 

Another kind of supersonic solution in the special case of 7 = 5/3 can be constructed 
assuming the outer boundary at some finite radius R out in which we require M. — 1. This 
solution takes the form 




(-) 



2/3 



(AH) 



where Rb = 2GM/cl ut and 5 = R out / Rb, and we denote p out and c out to be the density and 
sound speed at R ou t, respectively. Solution (All) exists for < 5 < 1/4. The dimensionless 
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accretion rate A = M / (AnR 2 B p out c OVi t) can be expressed in the form 

A = 5 2 . (A12) 

The Mach number in solution (All) 

^(^H^r <ai3) 

monotonically increases inward and approaches the asymptotic value Mo at R <C Rout- 
In the limit 5 <C 1, one gets Mo = 5~ 2 / 3 . In the case of marginal 5 = 1/4, one gets 
M. = Mo = 1 and solution (All) becomes equivalent to the self-similar solution (A9) in 
which a = l/y/2. Solution (All) approximates well our supersonic Model A4, assuming 
5 = 0.071 (short-dashed lines in Fig. 4). 



-22 - 



REFERENCES 

Abramowicz, M. A., Igumenshchev, I. V., Quataert, E., & Narayan, R. 2002, ApJ, 565, 1101 

Agol, E., & Kamionkowski, M. 2002, MNRAS, 334, 553 

Baganoff, F. K., et al. 2001, Nature, 413, 45 

Baganoff, F. K., et al. 2003, ApJ, 591, 891 

Belloni, T., Zampieri, L. & Campana, S. 1997, A&A, 319, 525 

Bisnovatyi-Kogan, G. S., & Ruzmaikin, A. A. 1974, Ap&SS, 28, 45 

Bisnovatyi-Kogan, G. S., & Lovelace, R. V. E. 1997, ApJ, 486, L43 

Bisnovatyi-Kogan, G. S., & Lovelace, R. V. E. 2000, ApJ, 529, 978 

Blaes, O., & Madau, P. 1993, ApJ, 403, 690 

Bondi, H. 1952, MNRAS, 112, 195 

Cabral, B., & Leedom, L. 1993, Computer Graphics: Proceedings: Annual Conference Series 
1993: SIGGRAPH 93 (New York: Association for Computing Machinery), 263 

Esin, A. A., Narayan, R., Ostriker, E., & Yi, I. 1996, ApJ, 465, 312 

Frank, J., King, A., & Raine, D. 1992, Accretion Power in Astrophysics (Cambridge Univ. 
Press, Cambridge) 

Fujita, Y., Inoue, S., Nakamura, T., Manmoto, T., & Nakamura, K. E. 1998, ApJ, 495, L85 

Ghez, A. M., et al. 2003, ApJ, 586L, 127 

Ho, L. C, Terashima, Y., & Okajima T. 2003, ApJ, 587L, 35 

Igumenshchev, I. V., & Abramowicz, M. A. 1999, MNRAS, 303, 309 

Igumenshchev, I. V., & Narayan, R. 2002, ApJ, 566, 137 (IN) 

Igumenshchev, I. V., Narayan, R., & Abramowicz, M. A. 2003, ApJ, 592, 1042 

Illarionov, A. F., & Sunyaev R. A. 1975, A&A, 39, 185 

Kato, S., Fukue, J., & Mineshige, S. 1998, Black-Hole Accretion Disks (Kyoto: Kyoto Univ. 
Press) 



-23- 



Katz J. L, 1977, ApJ, 215, 265 

Kovalenko, I. G., & Eremin, M. A. 1998, MNRAS, 298, 861 
Krumholz, M. R., McKee, C. F., & Klein, C. F. 2005, ApJ, 618, 757 

Landau, L. D., & Lifshitz, E. M. 1987, Electrodynamics of Continuous Media. Pergamon 
Press, Oxford 

Markovic, D. 1995, MNRAS, 277, 11 

Melia, F. 1992, ApJ, 387, L25 

Pen, U.-L., Matzner, C. D., & Wong, S. 2003, ApJ, 596, L207 

Perna, R., Narayan, R., Rybicki, G., Stella, L., & Treves, A. 2003, ApJ, 594, 936 

Narayan, R., Igumenshchev, I. V., & Abramowicz, M. A. 2000, ApJ, 539, 798 

Quataert, E. 2002, ApJ, 575, 855 

Quataert, E., & Gruzinov, A. 1999, ApJ, 520, 248 

Quataert, E., & Gruzinov, A. 2000, ApJ, 539, 809 

Shadmehri, M. 2004, A&A, 424, 379 

Shakura, N. I. 1972, AZh, 49, 921 

Shapiro, S. L., & Salpeter, E. E. 1975, ApJ, 198, 671 

Shvartsman, V. F. 1971, Soviet Astron. J., 15, 37 

Soria, R., Fabiano, G., Graham, A. W., Baldi, A., Elvis, M., Jerjen, H., Pellegrini, S., & 
Siemiginowska, A. 2006, ApJ, 40, 126 

Thorne, K. S., Price, R. H., & MacDonald, D. H. 1986, Black Holes: The Membrane 
Paradigm 

Toropin, Yu. M., Toropina, O. D., Savelyev, V. V., Romanova, M. M., Chechetkin, V. M., 
& Lovelace, R. V. E. 1999, ApJ, 517, 906 

Toropina, O. D., Romanova, M. M., Toropin, Yu. M., & Lovelace, R. V. E. 2003, ApJ, 953, 
472 

Treves, A., & Colpi, M. 1991, A&A, 241, 107 



-24- 

Turolla, R., Zampieri, L., Colpi, M., & Treves, A. 1994, ApJ, 426, L35 



This preprint was prepared with the AAS IATgX macros v5.2. 



-25 - 



Table 1. Parameters of the Models 



Model 


, a 

eo 




Reconnection heat 


Field topology 


Al 


1.2 









A2 


1.2 


10 2 




Toroidal 


A3 


1.2 


30 




Radial 


A4 


0.3 









Bl 


1.2 


10 


Included 


Loop 


B2 


1.2 


10 


Not included 


Loop 


CI 


0.3 


10 


Included 


Loop 


C2 


0.3 


10 


Not included 


Loop 



a Characterizes the specific internal energy of the injected mass given in 
units of e v ; r (see eq. [5]). The flows with e = 1.2 are subsonic and e = 0.3 
are supersonic. 

b Characterizes the initial or injected magnetic field strength and is de- 
fined in eq. (11). 
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Z A 




Fig. 1. — Geometry of the computational domain used in the simulations. The domain is 
extended from i? in to R out in the radial direction (R ou t/Rm = 10) and spans the angles 9 
and 0o in the polar and azimuthal directions, respectively (Oo = 0o = 7r/16). The source of 
gravity is located in the origin. 
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Injection of mass 

and magnetic field § 




boundary 

Fig. 2. — Schematic illustration of the simulation design. A three-dimensional MHD ac- 
cretion flow is formed as the result of a steady injection of mass and magnetic field in the 
vicinity of the outer radial boundary in the wedge computational domain. The geometry 
of the injected field is represented by magnetic lines (thick solid lines with arrows). The 
accretion flow with a frozen-in magnetic field is absorbed at the inner radial boundary. 
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Fig. 3. — Radial structure of the subsonic accretion flows in Models Al (short-dashed lines), 
Bl (solid lines), B2 (dotted lines), and self-similar solution (long-dashed lines, see eq. [A9]). 
Hydrodynamic Model Al is steady; turbulent MHD Models Bl and B2 are shown in quasi- 
steady states. All plotted quantities in the turbulent models — the density p, gas and mag- 
netic pressures P g and P m , temperature T, magnetic Mach number Ai m , accretion velocity 
v, and mass accretion rate M — have been averaged over the angles 9 and <fi and over the 
time interval r ~ 3£g. 
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Fig. 4. — Radial structure of the supersonic accretion flows in Models A4 (short-dashed 
lines), CI (solid lines), and C2 (dotted lines). Hydrodynamic Model A4 is steady; turbulent 
MHD Models Cf and C2 are shown in quasi-steady states. All plotted quantities in the 
turbulent models — the density p, gas and magnetic pressures P g and P m , temperature T, 
magnetic Mach number A4 m , accretion velocity v, and mass accretion rate M — have been 
averaged over the angles 9 and and over the time interval r ~ O.ltg. The MHD models 
have slowly moving outward Alfvenic shocks, which are located at the shown moment at 
R « 7.5-Rin in Model CI and R w 3i?i n in Model C2. 
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Fig. 5. — Snapshot of the density distribution in Model Bl in the equatorial plane. The 
source of gravity is located on the left, outside of the shown cross-section of the compu- 
tational domain (for a general view see Fig. 1). The shown radial range R out /Ri n = 10 
represents the fully simulated domain. The mass and magnetic field are steadily injected 
near the outer radial boundary on the right. The color bar on the left indicates the scale 
for logp (in arbitrary units). The model is turbulent and shown in a quasi-steady state. 
The fluctuations in the density are caused mainly by nonuniform reconnection heat and 
supported by convection motions. 

Fig. 6. — Snapshot of velocity streamlines in Model Bl in the equatorial plane in the same 
moment as in Fig. 5. The component of the streamlines parallel to the plane is shown. 
The flow pattern consists of the radially (horizontal in this view) extended narrow inflow- 
ing/outflowing streams and small-scale vortices. The streamlines have been plotted using 
the line integral convolution method of Cabral & Leedom (1993). 

Fig. 7. — Snapshot of magnetic lines in Model Bl in the equatorial plane in the same moment 
as in Fig. 5. The component of the field lines parallel to the plane is shown. One can see 
many magnetic loops stretched in the radial (horizontal in this view) direction. These loops 
are the result of convection and interchange instability. 

Fig. 8. — Snapshot of the distribution of the plasma (3 = PJ P m in Model Bl in the equatorial 
plane in the same moment as in Fig. 5. The color bar on the left indicates the scale for log/5. 
Small-scale regions of large (3, or weak magnetic fields, correspond to regions of reconnection 
and dissipation of magnetic energy. The large number of these regions indicates the high 
efficiency of the reconnection dissipation. Regions of low /3, or strong magnetic fields, are 
typically elongated in the radial (horizontal in this view) direction and associated with the 
inflowing streams. 

Fig. 9. — Snapshot of the distribution of specific entropy s in Model Bl in the equatorial 
plane in the same moment as in Fig. 5. The color bar on the left indicates the scale for 
logs (in arbitrary units). Regions of large s represent hot convective bubbles and streams, 
which have positive buoyancy and typically move outward. Regions of low s correspond to 
relatively cold inflowing matter. On average, the specific entropy is increased inward. 
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Fig. 10. — Snapshot of velocity streamlines in supersonic Model CI in the equatorial plane. 
The component of the streamlines parallel to the plane is shown. The flow pattern consists 
of two regions: the laminar super-fast-magnetosonic pre-shock (outer, on the right in this 
view) region and turbulent/convective post-shock (inner, on the left in this view) region. A 
quasi-steady Alfvenic shock separates these two regions and is located at R ~ 7 R ni in the 
shown moment. This shock slowly moves outward. 



Fig. 11. — Snapshot of magnetic lines in supersonic Model CI in the equatorial plane in the 
same moment as in Fig. 10. The component of the field lines parallel to the plane is shown. 
The field lines are purely radial (but oppositely directed in different sectors, see Fig. 2) and 
do not experience reconnections in the laminar pre-shock region (see capture to Fig. 10). 
The tangled field lines in the post-shock region are the result of convection and interchange 
instability. 



Fig. 12. — Radial distribution of the time-averaged total energy flux Fr normalized to the 
flux GMM- in / R in in Models B2 and C2 (solid lines). Reduction of Fr with respect to its 
value at the outer boundary (shown by dashed lines) represent losses of magnetic energy 
in reconnections. Supersonic Model C2 is shown at the moment, which corresponds to the 
Alfvenic shock position R ss 4.7 R m . Note the absence of energy losses in the laminar pre- 
shock region in Model C2. 
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